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segments. 



WO 03/060827 PCT/SG02/00006 



1 

METHOD AND APPARATUS FOR DETERMINING SYMMETRY 
IN 2D AND 3D IMAGES 

FIELD OF THE INVENTION 

The present invention relates to the processing of images. More 
5 particularly, the present invention relates to the processing of two-dimensional 
(2D) images or the processing of three-dimensional (3D) images, where such 3D 
images comprise a plurality of 2D image slices. In one aspect, the present 
invention is suitable for use in the determination of a plane of symmetry in a 3D 
medical image. For example, the present invention is suitable for use in 

10 determining the midsagittal plane (MSP) of the human brain from radiological 
images, and it will be convenient to hereinafter describe the invention in relation 
to that exemplary application. It should be appreciated, however, that the 
invention is not limited to that application, only. 
BACKGROUND OF THE INVENTION 

15 Human brains exhibit an approximate bilateral symmetry with respect to 

the inter-hemispheric (longitudinal) fissure bisecting the brain, known as the MSP. 
However, human brains are almost never perfectly symmetric (R. J. Davidson, K. 
Hugdahl, Eds., Brain Asymmetry. Cambridge, MA: MIT Press/Bradford Books, 
1996). A characteristic of the longitudinal fissure is that it is filled with cerebro- 

20 spinal fluid (CSF), which is quite dark in T1 -weighted magnetic resonance (MR) 
images and quite bright in T2-weighted MR images. 

MSP detection is often the first step in spatial normalization or anatomical 
standardization of brain images (J. L. Lancaster, T. G. Glass, B. R. Lankipalli, 
H. Downs, H. Mayberg, P. T. Fox, 'A modality-independent approach to spatial 

25 normalization of tomographic images of the human brain,' Human Brain Mapping 
(1995), 3: 209-223). It is also a useful first step in intrasubject inter/intramodality 
image registration (I. Kapouleas, A. Alavi, W. M. Alves, R. E. Gur, Registration of 
three-dimensional MR and PET images of the human brain without markers, 
Radiology (1991), 181 : 731-739). Determination or extraction of the MSP could 

30 provide a powerful tool to detect brain asymmetry due to tumors as well as any 
mass effects for diagnosis. In addition, extraction of the MSP is a prerequisite for 
automatically determining the anterior and posterior commissures needed by the 
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Talairach framework. The Talairach framework is a particular methodology that 
would be recognized by the person skilled in the art. 

Due to the great importance of the MSP, its determination and extraction 
has attracted quite a lot of research work. One paper entitled 'Robust midsagittal 
5 plane extraction from normal and pathological 3D neuroradiology images' IEEE 
Transactions on Medical Imaging (2001), 20(3): 175-192, by Y. Liu, R. T. Collins, 
W. E. Rothfus presented a cross-correlation algorithm to locate the plane that 
maximised bilateral symmetry. The algorithm first pre-processed the original 3D 
images to get edge maps. Then it calculated cross-correlation between axial 

10 slices and their reflection around an axis in the edge map to estimate the 
orientation that maximised the symmetry. From the estimated orientations of all 
the axial slices the orientation of the MSP was approximated. One of the 
limitations of the algorithm is its dependence on the edge map, which is quite 
prone to noise and skull appearances. Furthermore, the algorithm is hardly 

15 applicable to clinical environments since it is too time consuming (as long as 7 
minutes for a typical volume). The two papers by the same authors entitled 
'Evaluation of a robust midsagittal plane extraction algorithm for coarse, 
pathology 3D images' Proceedings of Medical Image Computing and Computer 
Assisted Intervention (MICCAI 2000): 81-94, and 'Automatic bilateral symmetry 

20 (midsagittal) plane extraction from pathology 3D neuroradiological images' 
Proceeding of SPIE International Symposium on Medical Imaging (1998), vol. 
3338: 1529-1539, respectively, have similar drawbacks. 

The paper entitled 'Automatic detection of mid-sagittal plane in 3D brain 
images' IEEE Transactions on Medical Imaging (1997), 16(6): 947-952 by B. A. 

25 Ardekani, J. Kershaw, M. Braun, I. Kanno presented a method purely based on 
symmetry assumption via calculation of cross correlation of grey level distribution 
around a plane, and took the plane with maximum cross-correlation as the 
approximated MSP. The algorithm was very sensitive to asymmetry in the 3D 
images, and could not deal with pathological images. The paper entitled 'Fully 

30 automatic identification of AC and PC landmarks on brain MRI using scene 
analysis' IEEE Transactions on Medical Imaging (1997), 16(5): 610-616 by L 
Verard, P. Allain, J. M. Travere, J. C. Baron, D. Bloyet presented a method to 
extract the MSP via detecting fissure lines in each axial slice. The algorithm 
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failed in those axial slices where the lateral ventricle was present. The algorithm 
also failed if any of the fissure lines were not straight or too broad, as in the case 
of excessive CSF. 

The paper entitled 'Hough transform detection of the longitudinal fissure in 
5 tomographic head images' IEEE Transactions on Medical Imaging (1991), 10(1): 
74-81 by M. E. Brummer presented a method to extract the MSP via Hough 
transform (HT) for plane detection. The algorithm needed pre-processing to get 
appropriate gradient modulus that should be fissure line dominant, which was 
difficult if not impossible due to noise, skull and skin appearances. The 

10 calculation was intensive for HT transformation in all slices. In addition, an 
incorrect detection of a fissure line in one slice would contribute to the final 
extraction, because there was no way to correct automatically. If the subjects 
being imaged were pathological such as containing excessive amount of CSF or 
having tumors, the algorithm also failed. 

15 Whilst the first method noted above by Liu et. al may deal with both normal 

and pathological images, none of the suggested algorithms are suitable for 
clinical application either due to excessively time intensive computing or due to an 
inability to process the real clinical data with pathology or ubiquitous asymmetry 
presented in axial slices. 

20 Essentially existing methods either try to extract the MSP via strict 

symmetry assumption without considering the grey level features of the fissure, or 
try to extract the MSP via fissure line detection based on low-level image 
processing techniques without considering the gross symmetry property of the 
fissure. None of the existing methods meet the current clinical requirements in 

25 terms of speed, accuracy, and robustness to clinical data. There is also a need to 
meet real time requirements for image processing in so far as determining the 
symmetry in an image. 
OBJECT OF THE INVENTION 

It is therefore an object of the present invention to provide a determination 

30 of symmetry in an image, which overcomes or ameliorates at least one of the 
problems of the prior art. 
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It is also an object of this invention to provide a clinically applicable 
solution to automatically determine the MSP of 3D radiological images. Ideally 
the solution should be specially tailored in terms of reasonably acceptable time 
for processing for a common hardware setup, robustness to variations in clinical 
5 data, user-friendliness, being fully automatic and stable. 
SUMMARY OF THE INVENTION 

With the above objects in mind, the present invention provides in one 
aspect a method of determining symmetry in an image, the method including the 
steps of: a) determine at least one searching line segment within a predefined 

10 search area of an image portion, the at least one searching line segment 
including a reference point (x, y) at its centre and an angle 6 with respect to a 
predetermined axis of the image portion; b) for each searching line segment, 
determine a first local characteristic in accordance with a measurement at points 
adjacent the searching line segment; c) determine the symmetry in the image in 

15 accordance with a calculation based on the first local characteristic. 

In essence the present invention stems from the realization that a localized 
characteristic measurement of a line segment in an image at points adjacent the 
line segment provides information at relatively low cost that can contribute to a 
relatively accurate determination of the symmetry of the entire image. A further 

20 realization of the present invention applies this localized characteristic 
measurement with histogram analysis to reject image anomalies that detract from 
a relatively accurate determination of the symmetry in an iterative manner that 
omits these anomalies from further processing, which would otherwise provide 
relatively low efficiency in image determination and analysis. In preferred 

25 embodiments, a further realization of the present invention combines the localized 
characteristic measurement with a further localised characteristic measurement of 
points along a given line segment in an image which is particularly useful in 3D 
image determination and processing. 



30 
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In one preferred form of the invention, the first local characteristic is a local 
symmetry index lsi(x,y,9) determined according to: d) for at least one pixel on said 
searching line segment, determining the sum of grey level differences for N pixel 
pairs located each side of said at least one pixel, said pixel pairs located on a line 
5 perpendicular to said searching line segment and intersecting said at least one 
pixel; e) summing said grey level differences determined in step d) and dividing 
this sum by the length of said searching line segment. 

In another preferred form of the invention, the first local characteristic is a 
modified local symmetry index mlsi(x,y,9) determined according to: f) for at least 

10 one pixel on said searching line segment, locating a closest pixel having a grey 
level less than a value below which a pixel is unlikely to be an object pixel and 
located on a line perpendicular to said searching line segment and passing 
through said pixel, and measuring a reference distance between said closest 
pixel and said line segment; g) for at least one pixel on said searching line 

15 segment, determining the sum of grey level differences for N pixel pairs located 
each side of said pixel, said pixel pairs located on said line perpendicular to said 
line segment and intersecting said at least one pixel; h) summing said grey level 
differences determined in step g) and dividing this sum by the length of said 
searching line segment 

20 In another aspect the present invention provides a method of determining 

a symmetry plane line segment of a 3D image, the 3D image comprising a set of 
parallel 2D slices, the method including the steps of: a) determine a first subset of 
slices from the plurality of 2D slices in accordance with predetermined criteria; b) 
for each of the first subset of slices, determine at least one searching line 

25 segment within a predefined search area of the slice, the at least one searching 
line segment including a reference point (x, y) at its centre and an angle 9 with 
respect to a predetermined axis of the slice; c) for each searching line segment, 
determine a first local characteristic in accordance with a measurement at points 
adjacent the searching line segment; d) for each searching line segment, 

30 determine a second local characteristic in accordance with a measurement at 
points along the searching line segment; e) determine the symmetry plane line 
segment in accordance with a calculation based on the first and second local 
characteristics. 
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In a further aspect the present invention provides a method of determining 
a symmetry line of a 2D image, the method including the steps of: a)determine at 
least one searching line segment within a predefined search area of a 2D image 
portion, the at least one searching line segment including a reference point (x, y) 
5 at its centre and an angle 9 with respect to a predetermined axis of the 2D image 
portion; b) for each searching line segment, determine a first local characteristic 
in accordance with a measurement at points adjacent the searching line segment; 
c) determine the symmetry line in accordance with a calculation based on the first 
local characteristic. 

10 In yet another aspect the present invention provides a method of 

calculating a symmetry plane in a 3D image, the 3D image comprising a plurality 
of parallel 2D slices, the method including the steps of: a) determine a first 
histogram of a first local characteristic of line segments within a predefined 
search area of selected slices, the first histogram including a peak first local 

15 characteristic; b) determine differences with respect to the peak first local 
characteristic for each first local characteristic; c) determine a first set of line 
segments having differences less than a first predetermined threshold; d) 
calculate an approximate symmetry plane in accordance with the first set of line 
segments; e) determine a second histogram of a second local characteristic of 

20 the first set of line segments, said second local characteristic being determined 
with reference to the approximate plane of symmetry, the second histogram 
including a peak second local characteristic; f) determine differences with respect 
to the peak second local characteristic for each second local characteristic; g) 
determine the end points of a second set of line segments having differences less 

25 than a second predetermined threshold; h) calculate the symmetry plane based 
on said end points of the second set of line segments. 

In another aspect the present invention provides a computer program 
product including: 

a computer usable medium having computer readable program code and 
30 computer readable system code embodied on said medium for determining 
symmetry in an image, said computer program product further including: 
computer readable code within said computer usable medium for: 



WO 03/060827 



PCT/SG02/00006 



7 

a) determining at least one searching line segment within a predefined 
search area of an image portion, the at least one searching line segment 
including a reference point (x, y) at its centre and an angle e with respect to a 
predetermined axis of the image portion; 
5 b) for each searching line segment, determining a first local 

characteristic in accordance with a measurement at points adjacent the searching 
line segment; 

c) determining the symmetry in the image in accordance with a 
calculation based on the first local characteristic. 
10 In a further aspect the present invention provides a computer program 

product including: 

a computer usable medium having computer readable program code and 
computer readable system code embodied on said medium for determining a 
symmetry plane line segment of a 3D image, the 3D image comprising a set of 
15 parallel 2D slices, said computer program product further including: 

computer readable code within said computer usable medium for: 

a) determining a first subset of slices from the plurality of parallel 2D 
slices in accordance with predetermined criteria; 

b) for each of the first subset of slices, determining at least one 
20 searching line segment within a predefined search area of the slice, the at least 

one searching line segment including a reference point (x, y) at its centre and an 
angle e with respect to a predetermined axis of the slice; 

c) for each searching line segment, determining a first local 
characteristic in accordance with a measurement at points adjacent the searching 

25 line segment; 

d) for each searching line segment, determining a second local 
characteristic in accordance with a measurement at points along the searching 
line segment; 

e) determining the plane line segment in accordance with a calculation 
30 based on the first and second local characteristics. 

In yet another aspect the present invention provides a computer program 
product including: 



WO 03/060827 



PCT/SG02/00006 



8 

a computer usable medium having computer readable program code and 
computer readable system code embodied on said medium for determining a 
symmetry line segment of a 2D image, said computer program product further 
including: 

5 computer readable code within said computer usable medium for: 

a) determining at least one searching line segment within a predefined 
search area of the 2D image, the at least one searching line segment including a 
reference point (x, y) at its centre and an angle 9 with respect to a predetermined 
axis of the slice; 

10 b) for each searching line segment, determining a first local 

characteristic in accordance with a measurement at points adjacent the searching 
line segment; 

c) determining the symmetry line segment in accordance with a 
calculation based on the first local characteristic. 
15 In another aspect the present invention provides a computer program product 
including: 

a computer usable medium having computer readable program code and 
computer readable system code embodied on said medium for calculating a 
symmetry plane in a 3D image, the 3D image comprising a plurality of parallel 2D 
20 slices, said computer program product further including: 

computer readable code within said computer usable medium for: 
a) determining a first histogram of a first local characteristic of line 
segments within a predefined search area of selected slices, the first histogram 
including a peak first local characteristic; 
25 b) determining differences with respect to the peak first local 

characteristic for each first local characteristic; 

c) determining a first set of line segments having differences less than 
a first predetermined threshold; 

d) calculating an approximate symmetry plane in accordance with the 
30 first set of line segments; 
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e) determining a second histogram of a second local characteristic of 
the first set of line segments, said second local characteristic being determined 
with reference to the approximate plane of symmetry, the second histogram 
including a peak second local characteristic; 
5 f) determining differences with respect to the peak second local 

characteristic for each second local characteristic; 

g) determining the end points of a second set of line segments having 
differences less than a second predetermined threshold; 

h) calculating the symmetry plane based on said end points of the 
10 second set of line segments. 

In another aspect the present invention provides a method of extracting 
MSP from radiological images including the steps of: determining axial slices to 
be processed from the radiological images; analysing the axial slices to determine 
fissure line segments; and calculating plane equation of MSP from the fissure line 
15 segments. 

The radiological images of the brain may be from different modalities (e.g. 
CT, MR, PET, SPET etc.,) or pulse sequences (e.g., T1 -weighted, T2-weighted, 
flair, or proton density weighted), arbitrary scanning orientations (axial, coronal, or 
sagittal), normal or abnormal brains. 

20 In the preferred embodiment, the fissure line segments may be determined 

by optimization of local symmetry index and fissure pattern measure. Further, 
histogram analysis may be used to calculate the plane equation. 

According to a further aspect, the present invention provides a computer 
program product including a computer usable medium having computer readable 

25 program code and computer readable system code embodied on said medium for 
extracting MSP from 3D radiological images within a data processing system, 
said computer program product further including computer readable code within 
said computer usable medium for: determining axial slices to be further 
processed for fissure line segments; extracting fissure line segments in axial 

30 slices by optimisation of local symmetry index and fissure pattern measure; and 
calculating plane equation of MSP from the extracted fissure line segments via 
histogram analysis. 
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In yet a further aspect, the present invention provides a method of 
calculating a plane equation from a plurality of non-intersecting line segments 
including the following steps: 

(a) calculating a first histogram of line angles of each line segment and 
5 finding a peak angle of the first histogram; 

(b) calculating angle differences between the peak angle and the line 
angles of each line segment and removing those line segments with the angle 
difference larger than a predetermined angle threshold; 

(c) approximating the plane equation from remaining line segments; 

10 (d) calculating distances of end points of the remaining line segments 

to the approximated plane and calculating a second histogram of distance 
distribution; 



(e) finding a peak distance of the second histogram; 

(f) calculating distance differences between the peak distance and the 



15 calculated distances, and removing those end points of the remaining line 
segments with the distance difference larger than a predetermined distance 
threshold; 

(g) calculating the plane equation based on the remaining end points of the 
remaining line segments. 
20 In still a further aspect, the present invention provides a method of 

determining a symmetry line segment in either a horizontal or vertical direction in 
a radiation image including the steps of: 



(d) calculating the symmetry line segment. 

Embodiments of the present invention are made possible using the novel 
concepts of local symmetry index, fissure pattern measure, and a novel way to 
calculate plane equation based on histogram analysis. The present invention 
30 makes use of the geometric property (approximate symmetry measured by local 
symmetry index) and may also make use of grey level (fissure pattern measure) 
properties of fissure line segments on axial slices for accurate and fast 
determination of fissure line segments. The robust approximation method to 



25 



(a) 
(b) 
(c) 



determining a searching region within the radiation image; 
calculating an approximate symmetry line segment; 
determining a refined searching region; 
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calculate the plane equation of MSP from extracted fissure line segments via 
histogram analysis allows the present invention to deal with pathological images 
correctly where either tumor(s) or mass effects dominate some of the axial slices, 
since the extracted fissure line segments in these abnormal slices will be 
5 regarded as outliers and are excluded from plane equation determination of MSP. 

According to another aspect, the present invention provides a method of 
extracting symmetry line from radiation images via maximizing the modified local 
symmetry index. 

According to yet a further aspect, the present invention provides a 
1 0 computer program product including a computer usable medium having computer 
readable program code and computer readable system code embodied on said 
medium for extracting symmetry line from radiation images within a data 
processing system, said computer program product further including computer 
readable code within said computer usable medium for extracting symmetry line 
1 5 from radiation images via maximising the modified local symmetry index. 
BRIEF DESCRIPTION OF THE DRAWINGS 

The invention may best be understood by reference to the following 
description in conjunction with the accompanying drawings, in which: 

Figure 1 shows a T1 -weighted MR slice with tumor, superimposed with the 
20 detected fissure line. 

Figure 2 shows the ideal head coordinate system X 0 Y 0 Z 0 and the real 
imaging coordinate system XYZ. 

Figure 3 is a flow chart illustrating the steps for extracting MSP from 3D 
radiological images according to an embodiment of the present invention. 
25 Figure 4 is a flow chart illustrating the steps for determining the axial slices 

to be processed further within given 3D radiological images according to an 
embodiment of the present invention. 

Figure 5 illustrates the concept of searching line segment ls(x,y,9). 
Figure 6 illustrates a pixel (x s ,y s ) on the searching line segment ls(x,y,0), 
30 and its associated 5 pairs of points on the opposite sides of the searching line 
segment, for calculating the contribution of the pixel (Xs,y s ) to local symmetry 
index. 
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Figures 7(a) and 7(b) show the fissure lines in a T1 -weighted MR slice and 
a T2-weighted MR slice. 

Figure 8 illustrates the steps to locate the fissure line segment on an axial 
slice according to an embodiment of the present invention. 
5 Figure 9 illustrates the steps to calculate the plane equation of MSP 

according to an embodiment of the present invention. 

Figure 10 illustrates how to form a closed polygon from all the remaining 
end points of fissure line segments. 

Figure 11 is an MSP calculated from a 3D T1 -weighted MR volume using 
10 the method disclosed in the present invention. 

Figure 12 is a flow chart illustrating steps for determining symmetry line 
from radiation images. 

Figure 13 shows a radiation image, superimposed with the detected 
symmetry line. 
15 DESCRIPTION OF PREFERRED EMBODIMENT 
Head Coordinate System 

Referring to Figure 2, the ideal head coordinate system is centered in the 
brain with positive X 0 , Y 0 , Z 0 axes pointing to the right, anterior, and superior 
directions, respectively. In clinical practice, the imaging coordinate system XYZ 
20 (Figure 2, black coordinate axis) differs from the ideal coordinate system due to 
positioning offsets (translations) and rotation of the head when imaged. With 
respect to this coordinate system, the plane equation X 0 = 0 is the MSP of the 
brain, and it is an objective of this invention to find the plane equation of MSP in 
the XYZ coordinate system. 
25 3D images and axial slices 

Neuroradiology scans are 3D volumetric data expressed as a stack of 2D 

slices. 

The 3D images can be obtained in 3 different ways: slices that are 
scanned along X, Y, and Z directions are called sagittal, coronal, and axial scans 
30 respectively. The scanned 3D images are linearly interpolated in space so that 
the voxel sizes are 1x1x1 mm 3 . Since from sagittal and coronal scans axial slices 
can be obtained, all the scans may be represented by axial slices. From here 
onwards, 3D images are isotropic in X, Y, and Z axes, with a voxel size being 
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1x1x1 mm 3 . The stack of 2D axial slices as a whole is also called a 3D volume. 
The 3D volume to be processed is denoted as Vol(x, y, z) with x, y, and z being 
the coordinates at voxel (x,y,z). Here x, y, and z are all non-negative integers 
satisfying 
5 0 < x < XSize 

0 < y < YSize 

0 < z < ZSize 

With these denotations, the number of axial slices is ZSize. Axial slices 
are counted from 0 to ZSize - 1 . The ith (0 < i < ZSize) axial slice is denoted as 

10 Si (So, Si Szsize-i), and it is the intersection of Vol(x,y,z) with the plane z = i. 

The grey level at a voxel (x,y,z) is denoted as g(x,y,z). In case of a 
radiation image, the grey level at pixel (x,y) is denoted as g(x,y). The two axes of 
an axial slice or a radiation image are denoted as H (horizontal) and V (vertical). 
Extracting MSP from axial slices 
15 One embodiment of the present invention will now be described in relation 

to extracting MSP from 3D MR images. 

For each axial slice Si, its intersection with MSP is a line that is physically 
called a fissure line. Since different axial slices are parallel, so the fissure lines 
are parallel as well. A fissure line segment is a line segment with two end points 
20 on the fissure line. 

Referring now to Figure 3, a method for extracting MSP once a scanned 
3D image has been read includes the following steps: 

1 . determining axial slices to be further processed for fissure line segments 
(10); 

25 2. approximating fissure line segments in axial slices by optimisation of local 
symmetry index and fissure pattern measure (20); and 
3. calculating the plane equation of MSP from the extracted fissure line 
segments via histogram analysis (60). 

The above-described steps which can provide MSP as an output will be 
30 described in greater detail herewith. 
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Determining the axial slices to be processed 

Now referring to Figure 4 for determining the axial slices to be processed 
further within a given 3D MR volume Vol(x, y, z). It is noted that some of the 
clinical data sets may have some slices vacant in the beginning or/and ending 

5 parts of the volume, which should be removed from further processing. Instead of 
having slices stretching from 0 to ZSize-1 , the first and last slices are determined 
by analysis of the histogram of the average grey levels (excluding starting/ending 
slices that have an average grey level smaller than a threshold described below). 
For each slice Si the system calculates the average grey level avg(Si) (11), then 

10 calculates the histogram of the average grey levels of all the axial slices (12). 
Once the histogram has been determined the system counts the accumulative 
percentage of the histogram from average grey level 0 to a threshold to, so that 
the accumulative percentage is equal to or above 30% (13). From axial slice 0 to 
slice ZSize -1 , the system finds the first slice with an average grey level equal to 

15 or bigger than the threshold t 0 and denotes the slice number as num 0 (14). 
Similarly, the system finds the last slice with an average grey level equal to or 
larger than the threshold t 0 and denotes the slice number as num n (15). From the 
beginning and ending slice numbers num 0 and num n , the slice increment is 
calculated to ensure the number of slices to be processed is smaller than 20 (16) 

20 as less than 20 slices has been found to provide a system which is not too time 
consuming but still maintains accuracy. The slice increment, which is denoted as 
StepSize, may be calculated as follows: 
L = num n - numo 
If (L > 30) StepSize = L/1 5 

25 Else if (L > 20) StepSize = L/1 0 

Else if (L > 15) StepSize = L/8 
Else StepSize = 1 . 

So the axial slices to be processed further will have the slice number snum 
calculated as follow: 
30 snum = num 0 + k * StepSize (k = 0, 1, L/StepSize) 
i.e., only S snu m will be processed further herewith. 
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Process axial slices to get the fissure line segments 

For the axial slice S snum to be processed further, its center of gravity 

(Xc(S S num), y c (S S num)) is calculated as follow: 

x c (S S num) = 22 x * g(x,y,snum)/ SZ g(x,y,snum) 
5 x y y 

yo(S snU m) = 22 y * g(x,y,snum)/ S£ g(x,y,snum) 

From a pixel (x,y) in slice S snU m, several line segments can be drawn. The 
line segment taking (x,y) as its center, with the length of line segment being 0.8 

10 YSize, and the angle with respect to the V axis being 9, is denoted as ls(x,y,e) as 
shown in Figure 5. Is(x,y,9) is called the searching line segment of pixel (x,y) with 
a searching angle 0, and pixel (x,y) is called the searching point. Note here that 
the rational to choose the length of line segment as 0.8 YSize is that the fissure 
line normally has a roughly vertical orientation, and its extension is around 0.8 

15 YSize, however the length of the line segment is not critical in the sense that a 
smaller length like 0.7 YSize can be used as well. All the searching points (x,y), 
together with the specified searching angle 6, satisfying the following conditions 
form the gross searching region: 

-10 mm < (x-Xc(Ssnum) ) ^ 10 mm 

20 y = yc(Ssnum) 

-1 0° < 9 < 1 0° 
where mm stands for millimeter. 

For a searching line segment ls(x,y,9), its local symmetry index lsi(x,y,9) is 
a measure of grey level symmetry in its local vicinity. Specifically, for each pixel 
25 (x s ,y s ) on the searching line segment, the system preferably checks 5 pairs of 
points at the opposite sides of ls(x,y,9) that are on the line perpendicular to 
ls(x,y,9) and passing (x,y) with a distance to ls(x,y,9) 2, 4, 6, 8 and 10 mm 
respectively (as shown in Figure 6). Denote cos(90°+9) as cos909, and 
sin(90°+9) as sin909. The contribution of the pixel (x s , y s ) to the local symmetry 
30 index lsi(x,y,0) is 

fabs ( g(x s +2cos909, y s +2stn909, snum) - g(x s -2cos909, y s -2sin909, 
snum) ) + fabs ( g(x s +4cos909, y s +4sin909, snum) - g(x s -4cos909, y s -4sin909, 
snum) ) + fabs ( g(x s +6cos909, y s +6sin909, snum) - g(x s -6cos909, y s -6sin909, 
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snum) ) + tabs ( g(x s +8cos909, y s +8sin909, snum) - g(x s -8cos909, y s -8sin909, 
snum) ) + tabs ( g(x s +10cos909, y s +10sin909, snum) - g(x s -10cos909, y s - 
10sin909, snum)) 

where fabs(.) is the absolute value function. For all (x s , y s ) on the searching line 
5 segment ls(x,y,9), the system calculates its contribution to lsi(x,y,9), then sums 
them up and divides the sum by the length of ls(x,y,9) to get the local symmetry 
index lsi(x,y,9), i.e., the formula to calculate lsi(x,y,9) is 
|ls(x,y,9)| lsi(x,y,9) = 

Jfajps ( g(x s +2cos909, y s +2sin909, snum) - g(x s -2cos909, y s -2sin909, snum) ) + 
E fabs ( g(x s +4cos909, y s +4sin909, snum) - g(x s -4cos909, y s -4sin909, snum) ) + 

(xs, ys) 

£ fabs ( g(x s +6cos909, y s +6sin909, snum) - g(x s -6cos909, y s -6sin909, snum) ) + 

(xs, ys) 

15 S fabs ( g(x s +8cos909, y s +8sin909, snum) - g(x s -8cos909, y s -8sin909, snum) ) + 

(xs, ys) 

Sfabs (g(x s +10cos909, y s +1 0sin909, snum) - g(x s -1 0cos909, y s -10sin909, snum)) 

(xs, ys) 

where |ls (x,y,9)| is the length of the searching line segment lsi(x,y,9). 

20 If the searching line segment is on the fissure line, then each pair of the 5 

pairs of points should have similar grey levels, hence the local symmetry index 
should be small. This means the local symmetry index introduced in this 
invention captures the grey level symmetry property in the vicinity of the 
searching line segment. Since the local vicinity of the searching line is around the 

25 center (from the definition of searching area) and is quite far from the skull and 
skin, it is insensitive to any asymmetry due to skull and/or skin. It is important to 
note that the asymmetry due to skull and skin is very frequent in clinical 
radiological images. Any inadequate handling of skull and skin could cause some 
deviation of the approximated fissure line from its correct fissure line as skull and 

30 skin in T1 -weighted MR slices will have bright grey level or strong edges and thus 
contribute to symmetry/asymmetry with a large weight if they are included in 
calculation for symmetry. None of the existing methods based on symmetry has 
excluded the influence of skin and skull. In addition, the implementation of local 
symmetry index is much faster in nature than the global symmetry index where all 
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the pairs of points on the opposite sides of the searching line segment with an 
equal distance are compared. 

Note that the local symmetry index does not consider the grey level feature 
of the fissure line itself. In fact, the fissure line also exhibits grey level features in 
5 normal slices. Figures 7(a) and 7(b) show the fissure lines in T1- and T2- 
weighted MR slices. In T1 -weighted MR slices, the fissure line is characterized 
by its low grey level, with most of the pixels having a lower grey level than the 
average grey level of the fissure line. Similarly, in T2-weighted slices, the fissure 
line is characterized by its large grey level, with most of the pixels having a higher 

10 grey level than the average grey level of the fissure line. In both cases, the 
standard deviation of the grey level on the fissure line is small. To quantify this 
kind of fissure, fissure pattern measure is introduced by the present invention. 
For each searching line segment ls(x,y,9), its fissure pattern measure is denoted 
as fpm(x,y,e). Denote the average grey level of the searching line segment 

15 ls(x,y,8) as agl(x,y,9), the standard deviation of grey level as sd(x,y,6), the 
number of pixels whose grey level is above the average as npa(x,y,6), and the 
number of pixels whose grey level is below average as npb(x,y,9), then the 
fissure pattern measure at the searching line segment ls(x,y,9) is calculated by: 
fpm(x,y,9) = 0.5(255-agl(x,y,9)) - sd(x,y,9) *0.8 YSize/ (3 npb(x,y,9) ) 

20 for T1 -weighted slices, and 

fpm(x,y,9) = 0.5agl(x,y,9) - sd(x,y,9) * 0.8 YSize/ (3 npa(x,y,9) ) 
for T2-weighted slices. 

Here the point is, if the searching line segment is on or close to the fissure 
line, the average grey level should be low in T1 -weighted MR slices and be high 

25 in T2-weighted MR slices, and the standard deviation of the grey level is low with 
most of the pixels being either above (T2-weighted) or below (T1 -weighted) the 
average grey level. As a result, the fissure pattern measure fpm(x,y,9) around the 
fissure line should be high. 

Referring to Figure 8, the steps to locate the fissure line includes: 

30 determination of gross searching region (21); calculation of local symmetry index 
lsi(x,y,9) for each searching line segment ls(x,y,9) in the gross searching region 
(22); calculation of fissure pattern measure fpm(x,y,9) (23); maximisation of the 
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weighted sum of local symmetry index and fissure pattern measure to determine 
the approximated fissure line segment (24). Here the weighted sum of local 
symmetry index and fissure pattern measure is denoted as ws(x,y,9), and is 
calculated by 
5 ws(x,y,9) = fpm(x,y,9) - W| S *lsi(x,y,9) 

where W| Si is a positive constant which could be decided by discriminant analysis 
and may vary in the range of 5.0 to 20.0. From here onwards, it is fixed as 10.0 
for the preferred system. 

The searching line segment ls(x,y,9) that maximises ws(x,y,9) is taken as 

10 the approximated fissure line segment at slice S sn um. The searching in the gross 
searching area is preferably carried out with x increment 2mm, and 9 increment 
2° in the first round to find the rough position of fissure line segment ls(x 0j yo,9 0 ). 
Then around the rough position of fissure line segment, the system defines the 
refined searching region such that for each (x,y,9) in the refined searching region, 

15 they satisfy 

-3.0 mm < (x-x 0 ) < 3.0 mm 

y = yo 

-3.0° < (9-9 0 ) <3.0° 

The search in the refined searching region is preferably done with x 
20 increment 0.5mm, and 9 increment 0.5°. The final extracted fissure line segment 
at slice S sn um is denoted by its two end points (x b (snum), y b (snum), snum), and 
(Xe(snum), y e (snum), snum) (where the first point is the beginning and the second 
point is the ending point of the line segment) and its angle with respect to V axis 
9(S S num). Since the extraction of fissure line segments may be wrong, all the 
25 extracted fissure line segments are called approximated fissure line segments 
herewith. 

Calculate the plane equation of MSP 

Referring to Figure 9, the calculation of the plane equation of MSP 
according to an embodiment is illustrated. 
30 The first step is to histogram the angles 9(S sn um) (61). Suppose the peak 

of the histogram is 9*. 
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The second step is to calculate the angle differences between 8(S snU m) and 
9*, and to remove those approximated fissure line segments that have an angle 
difference larger than 1.0 ° (62). The rational is, most of the approximated fissure 
line segments are correct and hence they should have very close orientation, 
5 while the approximated fissure lines in some of the abnormal slices will have a 
scattered orientation which should be excluded from calculation of plane equation 
of MSP. In this way, the outliers in orientation are removed. 

The third step is to approximate the plane equation from the remaining 
approximated fissure line segments (63). From those remaining approximated 
10 fissure line segments, the system forms a closed polygon by connecting 
beginning points of neighboring slices, connecting ending points of neighboring 
slices, and connecting the approximated fissure line segments in the two slices 
with smallest and largest z coordinates, as illustrated by Figure 10. The system 
may calculate the normal of the polygon by Newell's method (David Kirk, 
15 Graphics Gems III, p231-232, Academic Press, 1992). The approximated plane 
takes the calculated normal as the plane normal, and passes the average vector 
of all the remaining approximated fissure line segments. 

The fourth step is to calculate the distances of all the end points of the 
remaining approximated fissure line segments to the approximated plane and to 
20 create a histogram of the distance distribution (64). 

The fifth step is to find the peak of the histogram of the distance 
distribution (65). Denote the peak distance as d*. 

The sixth step is to calculate the distance differences between d* and all 
the calculated distances, and to remove those end points of the remaining 
25 approximated fissure line segments with a distance difference bigger than 1 .0 mm 
(66). After the removal based on orientation and distance histogram analysis, the 
remaining end points of the approximated fissure line segments are supposed to 
be the right end points of fissure line segments. 

The final step is to calculate the plane equation of MSP based on the 
30 remaining end points of the approximated fissure lines (67) by a similar method 
described in step 63: forming the closed polygon from remaining end points, 
approximating plane normal, and calculating the plane equation which takes the 
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calculated normal as the plane normal and passes the average vector of all the 
vertexes of the polygon. 

A working illustration of this embodiment of the invention is provided in 
Figure 1 1 . It shows the MSP of a T1 -weighted MR volume (256x256x1 68 voxels) 

5 extracted using the present embodiment of the invention. The correctness of the 
calculated MSP is assured via checking the extracted fissure line segments as 
well as via medical expert's knowledge of the expected MSP. The time taken for 
the extraction is 30 seconds on a SGI 02 system with mips 10000, 195 MHz 
CPU, and 190 MB RAM. 

10 Extract symmetry line in a radiation image 

The present invention can also be applied to extracting symmetry lines in 
radiation images. This is done via modified local symmetry index mlsi(x,y,9). 

Here the searching line segment ls(x,y,0) will have a segment length of 
0.8*XSize if the symmetry is approximately horizontal, or 0.8*YSize if the 

15 symmetry is approximately vertical. 

For a searching line segment ls(x,y,9), its modified local symmetry index is 
a measure of grey level symmetry in its local vicinity. Specifically, for each pixel 
(x s ,y s ) on the searching line segment, the system looks for the closest pixel which 
is on the line perpendicular to ls(x,y,8) and passing (x s ,y s ) (perpendicular line) 

20 with a grey level smaller than a value below which a pixel is unlikely to be an 
object pixel which could be decided by histogram analysis. Denote the distance 
between the closest pixel and pixel (Xs,y s ) as d 0 . Check 5 pairs of points on the 
perpendicular line with a distance to ls(x,y,8) d 0 +2, d 0 +4, d 0 +6, d 0 +8, and d 0 +10 
mm respectively. The contribution of the pixel (x s , y s ) to the modified local 

25 symmetry index mlsi(x,y,9) is 

fabs ( g(Xs+(do+2)cos909, y s +(d 0 +2)sin909) - g(x s -(d 0 +2)cos90e, 
y s (d 0 +2)sin909 ) ) + fabs ( g(x s +(d 0 +4)cos909, y s +(d 0 +4)sin909) - g(x s - 
(d 0 +4)cos909, y s -(d 0 +4)sin909 ) )+ fabs ( g(x s +(d 0 +6)cos909, y s +(d 0 +6)sin909) - 
g(x s -(d 0 +6)cos909, y s -(d 0 +6)sin909) )+ fabs ( g(x s +(d 0 +8)cos909, y s +(d 0 +8)sin909) 

30 - g(Xs-(d 0 +8)cos909, y s -(d 0 +8)sin909) ) + fabs ( g(x s +(d 0 +10)cos909, 
y s +(d 0 +10)sin909) - g(x s -(d 0 +1 0)cos909, y s -(d 0 +10)sin909) ) 
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For all (x s , y s ) on the searching line segment ls(x,y,9), calculate its 
contribution to mlsi(x,y,9), then sum them up and divide the sum by the length of 
ls(x,y,6) to get the modified local symmetry index lsi(x,y,9). i.e., the formula to 
calculate mlsi(x,y,G) is 
5 | ls(x,y,9)| mlsi(x,y,9) = 

E fabs (g(x s +(do+2)cos909, y s +(do+2)sin9O0)-g(x s -(d o +2)cos9O0, y s -(d o +2)sin9O0) ) + 

(xs, ys) 

2 fabs (g(x s +(d 0 +4)cos909, y s +(do+4)sin9Oe)-g(Xs-(do+4)cos9O0,y s -(do+4)sin9Oe) ) + 

(xs, ys) 

10 E fabs (g(x s +(do+6)cos9Oe,y s 4-(do+6)sin9O0)-g(x s -(do+6)cos9O0,y s -(do+6)sin9O0) ) + 

(xs, ys) 

E fabs(g(x s +(do+8)cos9O0,y s +(do+8)sin9O0)-g(x s -(do+8)cos9O0,ys-(do+8)sin9O0)) + 

(xs, ys) 

E fabs(g(x s +(d 0 +1 0)cos9O0,y s +(d o +1 0)sin9O0) - g(x s -(d 0 +1 0)cos9O0,y s -(d o +1 0)sin9O0)) 

^ g (xs, ys) 

Referring to Figure 12, the extraction of symmetry line from a radiation 
image according to an embodiment is illustrated. 

After reading in the radiation image, its gross searching region is 
determined (100). For approximately vertical symmetry, the gross searching 
20 region (x,y,9) is defined by 

-15 mm < (x-x c ) < 15 mm 

-2 mm < (y-y c ) < 2 mm 

-20° < 9 <20° 
where x c = XSize/2, and y c = YSize/2. 
25 For approximately horizontal symmetry, the gross searching region (x,y,9) 

is defined by: 

-2 mm < (x-Xc ) < 2 mm 

-15 mm < (y-y c ) < 15 mm 

70° < 9 < 110° 

30 The next step is to calculate the rough position of the symmetry line (101). 

This is achieved by minimizing the modified local symmetry index at searching 
points in the gross searching region. Search is done in the gross searching 
region preferably with x increment 2mm, y increment 2mm, and angle increment 
2°. The searching line segment ls(x,y,9) that minimizes the modified local 
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symmetry index is taken as the approximated symmetry line segment and is 

denoted as ls(x 0 ,yo,Oo). 

The next step is to determine the refined searching region (102). The 
refined region (x,y,9) is around (xo,yo,9o) such that: 
5 -3.0 mm < (x-x 0 ) < 3.0 mm 

-3.0 mm < (y-y 0 ) < 3.0 mm 
-3.0° < (6-e 0 ) < 3.0° 

The final step is to determine the final symmetry line (103). This is done 
via calculating the modified local symmetry index in the refined region and taking 

1 0 the line segment that minimizes the modified local symmetry index as the final 
symmetry line segment. The search in the refined searching region is preferably 
carried out with x increment 0.5mm, y increment 0.5mm, and 9 increment 0.5°. 

A working illustration of this embodiment of the invention is provided in 
Figure 13. It shows a thorax image with automatically computed symmetry axis 

15 superimposed. The time taken for the extraction is 50 seconds on a SGI 02 
system with mips 10000, 195 MHz CPU, and 190 MB RAM. 

This invention discloses a method for extracting MSP of human brain from 
radiological images. To gain speed, part of the axial slices in the volume are 
automatically picked out for fissure line segment extraction. Local symmetry 

20 index and fissure pattern measure are introduced to describe the grey level and 
geometrical features of fissure line segments in axial slices. The approximation 
of each fissure line segment is mathematically defined as locating the line 
segment that optimises the weighted sum of fissure pattern measure and local 
symmetry index. In order to tailor for clinical data and possible pathology, two ; 

25 rounds of outlier removal are applied. The first round is to remove those 
approximated fissure line segments that have quite a large angle deviation (1.0°) 
from the most probable angle of all the approximated fissure line segments. For 
the remaining approximated fissure line segments, a plane equation is 
approximated. The second round is to remove those end points of the remaining 

30 approximated fissure line segments that have quite a large distance deviation (1 .0 
mm) from the most probable distance to the approximated plane. Both the most 
probable angle and most probable distance are determined by finding the peak of 
the corresponding histograms. After two rounds of removal, the remaining end 
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points of the approximated fissure line segments are combined to calculate the 
plane equation of MSP. The disclosed method can extract MSP for clinical 
routine 3D MR images within 30 seconds under a common hardware setup 
regardless of the matrix dimension of the images (as large as 256x256x256), 
5 pulse sequences (T1 -weighted, T2-weighted, flair, or proton density weighted), 
scanning orientations (axial, coronal or sagittal), and normal or pathological 
subjects. 

The preferred embodiment includes 3 major steps: (1) determining axial 
slices to be further processed for fissure line segments, (2) approximating fissure 
1 o line segments in axial slices by optimisation of local symmetry index and fissure 
pattern measure, and (3) calculating the plane equation of MSP from the 
approximated fissure line segments via histogram analysis. 

With slight modification, the method is applied to extracting symmetry line 
in a radiation image. 

15 As the present invention may be embodied in several forms without 

departing from the spirit of the essential characteristics of the invention, it should 
be understood that the above described embodiments are not to limit the present 
invention unless otherwise specified, but rather should be construed broadly 
within the spirit and scope of the invention as defined in the appended claims. 

20 Various modifications and equivalent arrangements are intended to be included 
within the spirit and scope of the invention and appended claims. 
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CLAIMS: 

1 . A method of determining symmetry in an image, the method including the 
steps of: 

a) determine at least one searching line segment within a predefined 
search area of an image portion, the at least one searching line segment 
including a reference point (x, y) at its centre and an angle 9 with respect to a 
predetermined axis of the image portion; 

b) for each searching line segment, determine a first local 
characteristic in accordance with a measurement at points adjacent the searching 
line segment; 

c) determine the symmetry in the image in accordance with a 
calculation based on the first local characteristic. 

2. A method as claimed in claim 1 wherein, the first local characteristic is a 
local symmetry index lsi(x,y,9) determined according to: 

a) for at least one pixel on said searching line segment, determining 
the sum of grey level differences for N pixel pairs located each side of said at 
least one pixel, said pixel pairs located on a line perpendicular to said searching 
line segment and intersecting said at least one pixel; 

b) summing said grey level differences determined in step d) and 
dividing this sum by the length of said searching line segment. 

3. A method as claimed in claim 1 wherein, the first local characteristic is a 
modified local symmetry index mlsi(x,y,9) determined according to: 

a) for at least one pixel on said searching line segment, locating a 
closest pixel having a grey level less than a value below which a pixel is unlikely 
to be an object pixel and located on a line perpendicular to said searching line 
segment and passing through said pixel, and measuring a reference distance 
between said closest pixel and said line segment; 
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b) for at least one pixel on said searching line segment, determining 
the sum of grey level differences for N pixel pairs located each side of said pixel, 
said pixel pairs located on said line perpendicular to said line segment and 
intersecting said at least one pixel; 

c) summing said grey level differences determined in step g) and 
dividing this sum by the length of said searching line segment. 

4. A method as claimed in claim 2 wherein, N = 5 and the pixel pairs are 
located at 2, 4, 6, 8 and 1 0mm, respectively, from said at least one pixel. 

5. A method as claimed in claim 3 wherein, N = 5 and the pixel pairs are located 
at 2, 4, 6, 8 and 10mm plus said reference distance, respectively, from said at 
least one pixel. 

6. A method of determining a symmetry plane line segment of a 3D image, 
the 3D image comprising a set of parallel 2D slices, the method including the 
steps of: 

a) determine a first subset of slices from the plurality of 2D slices in 
accordance with predetermined criteria; 

b) for each of the first subset of slices, determine at least one 
searching line segment within a predefined search area of the slice, the at least 
one searching line segment including a reference point (x, y) at its centre and an 
angle e with respect to a predetermined axis of the slice; 

c) for each searching line segment, determine a first local 
characteristic in accordance with a measurement at points adjacent the searching 
line segment; 

d) for each searching line segment, determine a second local 
characteristic in accordance with a measurement at points along the searching 
line segment; 

e) determine the symmetry plane line segment in accordance with a 
calculation based on the first and second local characteristics. 



WO 03/060827 



PCT/SG02/00006 



26 

7. A method as claimed in claim 6 wherein, the first local characteristic is a 
local symmetry index lsi(x,y,0) determined according to: 

a) for at least one pixel on said searching line segment, determining 
the sum of grey level differences for N pixel pairs located each side of said at 
least one pixel, said pixel pairs located on a line perpendicular to said searching 
line segment and intersecting said at least one pixel; 

b) summing said grey level differences determined in step f) and 
dividing this sum by the length of said searching line segment. 

8. A method as claimed in claim 7 wherein, N = 5 and the pixel pairs are 
located at 2, 4, 6, 8 and 10mm, respectively, from said at least one pixel. 

9. A method as claimed in claim 6, 7 or 8 wherein, the second local 
characteristic is a fissure pattern measure fpm(x,y,0) determined according to: 

fpm(x,y,9) = 0.5(255-agl(x,y, 0)) - sd(x,y,9) *0.8 YSize/ (3 npb(x,y,0) ) 
for T1 -weighted slices, and 

fpm(x,y,9) = 0.5agl(x,y, 9) - sd(x,y,9) * 0.8 YSize/ (3 npa(x,y,9) ) 
for T2-weighted slices, 

where agl(x,y,9) is the average grey level along said searching line 
segment ls(x,y,9), sd(x,y,9) the standard deviation of grey level along Is(x,y,9), 
npa(x,y,9) the number of pixels on ls(x,y,9) whose grey level is above agi(x,y,9), 
and npb(x,y,9) the number of pixels on ls(x,y,9) whose grey level is below 
agl(x,y). 

10. A method as claimed in claim 6, 7, 8 or 9 wherein, the calculation of step 
e) includes a weighted sum ws(x,y,9) of the first and second local characteristics. 

11. A method as claimed in claim 10 wherein, said weighted sum ws(x,y,0) of 
said local symmetry index lsi(x,y,9) and said fissure pattern measure fpm(x,y,9) 
along said searching line segment lsi(x,y) is defined by: 

ws(x,y,0) = ai* fpm(x,y,9) - a 2 * lsi(x,y,9) 
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where and a 2 are positive floating constants that are determined by 
discriminant analysis. 

12. A method as claimed in any one of claims 6 to 11 wherein the method is 
performed in two steps of different resolutions including: 

(a) searching in a first predefined gross search area of a slice defined 

by 

-10 mm < (x-Xc(Ssnum) ) ^ 10 mm 

y ~ yc(Ssnum) 

-10° < 0 < 10° 

with a gross searching increment, given by an x increment of 2mm, and an 
angle increment of 2°, to estimate the position of said symmetry plane line 
segment, where x c (S sn um) and y c (S S num) represent the centre of said slice; 

(b) searching in a second predefined fine search area of a slice defined 

by 

-3.0 mm < (x-x 0 ) < 3.0 mm 
y = Yo 

-3.0° < (9-e 0 ) < 3.0° 
with a fine searching increment, given by an x increment of 0.5mm, and an 
angle increment of 0.5°, where (x 0 , y 0 , e 0 ) represents the searching line segment 
that maximises said weighted sum of said local symmetry index lsi(x,y,0) and 
said fissure pattern measure fpm(x,y,9) in said gross search area with said gross 
searching increment. 

13. A method as claimed in any one of claims 6 to 12 wherein said first subset 
of slices from the set of 2D slices are determined in accordance with the following 
predetermined criteria: 

(a) calculate average grey levels for each slice; 

(b) calculate a histogram of said average grey levels; 

(c) determine a grey level threshold; 

(d) determine a first slice with a grey level equal to or greater than said 
threshold; 
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(e) determine a last slice with a grey level equal to or greater than said 
threshold; and 

(f) determine a slice increment step size between said first slice and 
said last slice. 

14. A method as claimed in any one of claims 6 to 13 wherein the symmetry 
plane line segment is a fissure line segment of a medical image. 

15. A method as claimed in claim 14 wherein, the medical image is a 
volumetric image produced by any one or more of the following techniques: 

a) CT; 

b) MRI; 

c) PET; 

d) SPET. 

16. A method of determining a symmetry line of a 2D image, the method 
including the steps of 

a) determine at least one searching line segment within a predefined 
search area of a 2D image portion, the at least one searching line segment 
including a reference point (x, y) at its centre and an angle G with respect to a 
predetermined axis of the 2D image portion; 

b) for each searching line segment, determine a first local 
characteristic in accordance with a measurement at points adjacent the searching 
line segment; 

c) determine the symmetry line in accordance with a calculation based 
on the first local characteristic. 

17. A method as claimed in claim 16 wherein, the first local characteristic is a 
modified local symmetry index mlsi(x,y,9) determined according to: 

(a) for at least one pixel on said searching line segment, locating a 
closest pixel having a grey level less than a value below which a pixel is unlikely 
to be an object pixel and located on a line perpendicular to said searching line 
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segment and passing through said pixel, and measuring a reference distance 
between said closest pixel and said line segment; 

(b) for at least one pixel on said searching line segment, determining 
the sum of grey level differences for pixel pairs located each side of said pixel, 
said pixel pairs located on said line perpendicular to said line segment and 
intersecting said at least one pixel; 

c) summing said grey level differences determined in step b) and 
dividing this sum by the length of said searching line segment. 

18. A method as claimed in any one of claims 16 or 17 wherein the symmetry 
line is a symmetry line of a medical image. 

19. A method as claimed in claim 18 wherein, the medical image is a 2D 
radiation image. 

20. A method of calculating a symmetry plane in a 3D image, the 3D image 
comprising a plurality of parallel 2D slices, the method including the steps of: 

a) determine a first histogram of a first local characteristic of line 
segments within a predefined search area of selected slices, the first histogram 
including a peak first local characteristic; 

b) determine differences with respect to the peak first local 
characteristic for each first local characteristic; 

c) determine a first set of line segments having differences less than a 
first predetermined threshold; 

d) calculate an approximate symmetry plane in accordance with the 
first set of line segments; 

e) determine a second histogram of a second local characteristic of the 
first set of line segments, said second local characteristic being determined with 
reference to the approximate plane of symmetry, the second histogram including 
a peak second local characteristic; 

f) determine differences with respect to the peak second local 
characteristic for each second local characteristic; 



WO 03/060827 



PCT/SG02/00006 



30 

g) determine the end points of a second set of line segments having 
differences less than a second predetermined threshold; 

h) calculate the symmetry plane based on said end points of the 
second set of line segments. 

21 . A method as claimed in claim 20 wherein, the first local characteristic is the 
angle of a line segment with respect to a predetermined axis of the selected 
slices and the first predetermined threshold is an angle threshold of 1 .0°. 

22. A method as claimed in claim 20 or 21 wherein, the second local 
characteristic is the distance from end points of line segments of the first set of 
line segments to said approximate symmetry plane and the second 
predetermined threshold is a distance threshold of 1.0 unit length. 

23. A method as claimed in claim 20, 21 or 22 wherein, the line segments are 
symmetry plane line segments as claimed in claims 6 to 19. 

24. A method as claimed in any one of claims 20 to 23 wherein the symmetry 
plane is a midsagittal plane MSP of the human brain. 

25. A computer program product including: 

a computer usable medium having computer readable program code and 
computer readable system code embodied on said medium for determining 
symmetry in an image, said computer program product further including: 

computer readable code within said computer usable medium for: 

a) determining at least one searching line segment within a predefined 
search area of an image portion, the at least one searching line segment 
including a reference point (x, y) at its centre and an angle 9 with respect to a 
predetermined axis of the image portion; 

b) for each searching line segment, determining a first local 
characteristic in accordance with a measurement at points adjacent the searching 
line segment; 
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c) determining the symmetry in the image in accordance with a 
calculation based on the first local characteristic. 

26. A computer program product including: 

a computer usable medium having computer readable program code and 
computer readable system code embodied on said medium for determining a 
symmetry plane line segment of a 3D image, the 3D image comprising a set of 
parallel 2D slices, said computer program product further including: 

computer readable code within said computer usable medium for: 

a) determining a first subset of slices from the plurality of 2D slices in 
accordance with predetermined criteria; 

b) for each of the first subset of slices, determining at least one 
searching line segment within a predefined search area of the slice, the at least 
one searching line segment including a reference point (x, y) at its centre and an 
angle 9 with respect to a predetermined axis of the slice; 

c) for each searching line segment, determining a first local 
characteristic in accordance with a measurement at points adjacent the searching 
line segment; 

d) for each searching line segment, determining a second local 
characteristic in accordance with a measurement at points along the searching 
line segment; 

e) determining the symmetry plane line segment in accordance with a 
calculation based on the first and second local characteristics. 

27. A computer program product including: 

a computer usable medium having computer readable program code and 
computer readable system code embodied on said medium for determining a 
symmetry line segment of a 2D image, said computer program product further 
including: 

computer readable code within said computer usable medium for: 
a) determining at least one searching line segment within a predefined 
search area of the 2D image, the at least one searching line segment including a 
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reference point (x, y) at its centre and an angle 9 with respect to a predetermined 
axis of the slice; 

b) for each searching line segment, determining a first local 
characteristic in accordance with a measurement at points adjacent the searching 
line segment; 

c) determining the line segment in accordance with a calculation 
based on the first local characteristic. 

28. A computer program product including: 

a computer usable medium having computer readable program code and 
computer readable system code embodied on said medium for calculating a 
symmetry plane in a 3D image, the 3D image comprising a plurality of parallel 2D 
slices, said computer program product further including: 

computer readable code within said computer usable medium for: 

a) determining a first histogram of a first local characteristic of line 
segments within a predefined search area of selected slices, the first histogram 
including a peak first local characteristic; 

b) determining differences with respect to the peak first local 
characteristic for each first local characteristic; 

c) determining a first set of line segments having differences less than 
a first predetermined threshold; 

d) calculating an approximate symmetry plane in accordance with the 
first set of line segments; 

e) determining a second histogram of a second local characteristic of 
the first set of line segments, said second local characteristic being determined 
with reference to the approximate plane of symmetry, the second histogram 
including a peak second local characteristic; 

f) determining differences with respect to the peak second local 
characteristic for each second local characteristic; 

g) determining the end points of a second set of line segments having 
differences less than a second predetermined threshold; 

h) calculating the symmetry plane based on said end points of the 
second set of line segments. 
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29. A method of extracting MSP of a brain from radiological images including 
the steps of: 

(a) determining axial slices to be processed from said radiological 
images; 

(b) analysing said axial slices to determine fissure line segments; and 

(c) calculating plane equation of MSP from said fissure line segments. 

30. A method as claimed in claim 29 wherein the said radiological images of 
brain can be from different modalities (e.g. CT, MR, PET, SPET etc.,) or pulse 
sequences (e.g., T1 -weighted, T2-weighted, flair, or proton density weighted), 
arbitrary scanning orientations (axial, coronal, or sagittal), normal or abnormal 
brains. 



31. A method as claimed in claim 29 or claim 30, wherein said fissure line 
segments are determined by optimisation of local symmetry index and fissure 
pattern measure. 

32. A method as claimed in any one of claims 29 to 31 , wherein histogram 
analysis is used to calculate said plane equation. 



33. A method as claimed in any one of claims 29 to 32 wherein said axial 
slices to be processed are determined using the following steps: 

(a) calculating average grey levels for each slice; 

(b) calculating a histogram of said average grey levels; 

(c) determining grey level threshold; 

(d) determining a first slice with grey level equal to or greater than said 
threshold; 

(e) determining a last slice with grey level equal to or greater than said 
threshold; and 

(f) determining slice increment step size between said first slice and 
said last slice. 
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34. A method as claimed in claim 33 wherein said threshold is equal to the 
smallest grey level with an accumulative percentage of said histogram equal to or 
greater than 30%. 

35. A method as claimed in claim 33 or claim 34, wherein said slice increment 
step size is determined such that 20 or less axial slices are processed. 

36. A method as claimed in any one of claims 29 to 35 further including the 
process of determining said fissure line segments via the steps of: 

(a) calculating local symmetry index along searching line segments; 

(b) calculating fissure pattern measure along said searching line 
segments; 

(c) locating said fissure line segment by finding said searching line 
segment with maximum weighted sum of said local symmetry index and said 
fissure line pattern. 

37. A method as claimed in claim 36, wherein said local symmetry is 
calculated by: 

(a) for any pixel on said searching line segment, summing grey level 
differences for pixel pairs located either side of said pixel, said pixel pairs being 
on a line perpendicular to said line segment and passing through said pixel; 

(b) summing grey level differences of each said pixel on said line 
segment and dividing this sum by the length of said searching line segment. 

38. A method as claimed in claim 37, wherein five said pixel pairs are 
summed. 

39. A method as claimed in claim 38 wherein said pixel pairs are located 2, 4, 
6, 8 and 10 mm from said pixel. 

40. A method as claimed in claim 36, wherein said fissure pattern measure 
fpm(x,y,0) along said searching line segment ls(x,y,0) is calculated by: 

fpm(x,y,0) = 0.5(255-agI(x,y, 0)) - sd(x,y,0) *0.8 YSize/ (3 npb(x,y,0) ) 
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for T1 -weighted slices, and 

fpm(x,y,0) = 0.5agl(x,y, G) - sd(x,y,0) * 0.8 YSize/ (3 npa(x,y,0) ) 
for T2-weighted slices, 

where agl(x,y,8) is the average grey level along said line segment ls(x,y,0), 
sd(x,y,0) the standard deviation of grey level along ls(x,y,9), npa(x,y,0) the 
number of pixels on ls(x,y,0) whose grey level is above agl(x,y,0), and npb(x,y,0) 
the number of pixels on ls(x,y,0) whose grey level is below agl(x,y). 

41 . A method as claimed in any one of claims 36 to 40, wherein said weighted 
sum ws(x,y,0) of said local symmetry index lsi(x,y,0) and said fissure pattern 
measure fpm(x,y,0) along said searching line segment lsi(x,y) is defined by: 

ws(x,y,0) = a!* fpm(x,y,0) - a 2 * lsi(x,y,0) 

where ai and a 2 are positive floating constants that are determined by 
discriminant analysis. 

42. A method as claimed in any one of claims 36 to 41 , wherein determining 
said fissure line segments is completed in two steps of different resolutions, 
namely 

(a) searching in a broad region defined by 
-10 mm < (x-Xc(Ssnum) ) < 10 mm 

y = yc(Ssnum) 

-10° <0 <10° 

with gross searching increment (x increment 2mm, and angle increment 
2°) to estimate the position of said fissure line segment, where x c (S sn um) and 
yc(Ssnum) represent the centre of said slice; 

(b) searching in a refined region defined by 
-3.0 mm < (x-x 0 ) < 3.0 mm 

y = yo 

-3.0° < (0-0 o ) < 3.0° 
with x increment 0.5mm, and angle increment 0.5°, where (x 0 , yo, ©o) 
represents said searching line segment that maximises said weighted sum of 
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said local symmetry index and said fissure pattern measure in said broad region 
with said gross searching increment. 

43. A method as claimed in any one of claims 29 to 42 wherein said plane 
equation of MSP is calculated from approximated fissure line segments via 
histogram analysis including the steps of: 

(a) calculating an angle histogram of angles of fissure line segments for 
said axial slices and determining peak angle of said angle histogram; 

(b) calculating angle differences between said angles and said peak 
angle, and removing those said approximated fissure line segments having an 
angle difference larger than 1 .0°; 

(c) approximating plane equation from the remaining said approximated 
fissure line segments; 

(d) calculating distances of all end points of said remaining 
approximated fissure line segments to said approximated plane and calculating 
the distance histogram of the distance distribution; 

(e) finding the peak distance of said distance histogram of distance 
distribution; 

(f) calculating the distance differences between said peak distance and 
all said calculated distances, and removing those end points of said approximated 
fissure line segments with a distance difference larger than 1.0 mm; 

(g) calculating said plane equation of MSP based on remaining end 
points of said approximated fissure line segments. 

44. A method of calculating a plane equation from a plurality of non- 
intersecting line segments including the following steps: 

(a) calculating a first histogram of line angles of each said line segment 
and finding a peak angle of said first histogram; 

(b) calculating angle differences between said peak angle and said line 
angles of each said line segment and removing those said line segments with 
said angle difference larger than a predetermined angle threshold; 

(c) approximating said plane equation from remaining line segments; 
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(d) calculating distances of end points of said remaining line segments 
to said approximated plane and calculating a second histogram of distance 
distribution; 

(e) finding a peak distance of said second histogram; 

(f) calculating distance differences between said peak distance and 
said calculated distances, and removing those end points of said remaining line 
segments with said distance difference larger than a predetermined distance 
threshold; 

(g) calculating said plane equation based on the remaining end points 
of said line segments. 

45. A method as claimed in claim 44, wherein said angle threshold is 1 .0°. 

46. A method as claimed in claim 44 or claim 45, wherein said distance 
threshold is 1 .0 unit length. 

47. A method of determining a final symmetry line segment in either a 
horizontal or vertical direction in a radiation image including the steps of: 

(a) determining a gross searching region within said radiation image; 

(b) calculating an approximate symmetry line segment; 

(c) determining a refined searching region; 

(d) calculating said final symmetry line segment. 

48. A method as claimed in claim 47, wherein said gross searching region is 
determined by: 

for vertical symmetry 
-15 mm < (x-Xc ) < 15 mm 
-2 mm < (y-y c ) < 2 mm 
-20° <0 <20° 
for horizontal symmetry 
-2 mm < (x-Xc ) < 2 mm 
-15mm<(y-y c )<15mm 
70° <0 <110° 
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where Xc = XSize/2, and y c = YSize/2, and XSize and YSize are the 
maximum X and Y coordinates. 

49. A method as claimed in claim 47 or claim 48 wherein said approximate 
symmetry line segment is determined by minimizing modified local symmetry 
index in said gross searching region with an x increment 2 mm, y increment 2 
mm, and 9 increment 2°. 

50. A method as claimed in claim 49 wherein said modified local symmetry 
index of a searching line segment is defined by: 

(a) for each pixel on said searching line segment, locating a closest 
pixel having a grey level less than a value below which a pixel is unlikely to be an 
object pixel which could be decided by histogram analysis and located on a line 
perpendicular to said searching line segment and passing through said pixel, and 
measuring a reference distance between said closest pixel and said line segment; 

(b) for each said pixel on said searching line segment, summing grey 
level differences for pixel pairs located either side of said pixel, said pixel pairs 
being on said line perpendicular to said line segment and passing through said 
pixel; 

(c) summing grey level differences of each said pixel on said line 
segment and dividing this sum by the length of said searching line segment. 

51. A method as claimed in claim 50, wherein five said pixel pairs are 
summed. 

52. A method as claimed in claim 51 wherein said pixel pairs are located said 
reference distance plus 2, 4, 6, 8 and 1 0mm respectively from said pixel. 

53. A method as claimed in any one of claims 47 to 52, wherein said refined 
searching region is defined by 

-3.0 mm < (x-x 0 ) < 3.0 mm 
-3.0 mm < (y-y 0 ) < 3.0 mm 
-3.0° < (9-e 0 ) < 3.0° 
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54. A method as claimed in any one of claims 50 to 53 wherein said final 
symmetry line segment is determined by the steps of: 

(a) calculating said modified local symmetry index in said refined 
searching region; 

(b) searching in said refined region with x increment 0.5mm, y 
increment 0.5mm, and 9 increment 0.5°; 

(c) locating said searching line segment with minimised said modified 
local symmetry index. 

55. A computer program product including a computer usable medium having 
computer readable program code and computer readable system code embodied 
on said medium for extracting MSP from 3D radiological images within a data 
processing system, said computer program product further including computer 
readable code within said computer usable medium for: determining axial slices 
to be processed for fissure line segments; extracting fissure line segments in axial 
slices by optimisation of local symmetry index and fissure pattern measure; and 
calculating plane equation of MSP from said extracted fissure line segments via 
histogram analysis. 

56. A product as claimed in claim 55, wherein said axial slices to be processed 
are determined using the following steps: 

(a) calculating average grey levels for each slice; 

(b) calculating a histogram of said average grey levels; 

(c) determining a grey level threshold which is the smallest grey level 
with an accumulative percentage of said histogram equal to or greater than 30%; 

(d) determining a first slice with grey level equal to or greater than said 
threshold; 

(e) determining a last slice with grey level equal to or greater than said 
threshold; and 

(f) determining slice increment step size between said first slice and 
said last slice, such that 20 or less axial slices are processed. 
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57. A product as claimed in any one of claims 55 to 56 further including the 
process of determining said fissure line segments via the steps of: 

(a) calculating local symmetry index along searching line segments; by 

(i) for any pixel on said searching line segment, summing grey 
level differences for five pixel pairs located 2, 4, 6, 8 and 10 mm from said 
pixel and on either side of said pixel, said pixel pairs being on a line 
perpendicular to said line segment and passing through said pixel; 

(ii) summing grey level differences of each said pixel on said line 
segment and dividing this sum by the length of said searching line 
segment; 

(b) calculating fissure pattern measure along said searching line 
segments; 

(c) locating said fissure line segment by finding said searching line 
segment with maximum weighted sum of said local symmetry index and said 
fissure line pattern. 

58. A method of extracting symmetry line from radiation images via maximising 
the modified local symmetry index. 
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